###############################################
#FIGURE C1
###############################################

library(estimatr)
library(ggplot2)
library(grid)

####load cleaned data####

data1 <- read.csv("./data/data1.csv")
data2 <- read.csv("./data/data2.csv")
data3 <- read.csv("./data/data3.csv")


####results without controls####

lm1 <- lm_robust(list_items ~ as.factor(trt)*sens, data=data1)
lm2 <- lm_robust(list_items ~ as.factor(trt)*sens, data=data2)
lm3 <- lm_robust(list_items ~ as.factor(trt)*sens, data=data3)

com1 <- c("Identity Threat\nvs. Control", lm1$coefficients[4], (lm1$coefficients[4] + (1.65*lm1$std.error[4])), (lm1$coefficients[4] - (1.65*lm1$std.error[4])),
          (lm1$coefficients[4] + (1.96*lm1$std.error[4])), (lm1$coefficients[4] - (1.96*lm1$std.error[4]))) 
com2 <- c("Whataboutism\nvs. Identity Threat", lm2$coefficients[4], (lm2$coefficients[4] + (1.65*lm2$std.error[4])), (lm2$coefficients[4] - (1.65*lm2$std.error[4])),
          (lm2$coefficients[4] + (1.96*lm2$std.error[4])), (lm2$coefficients[4] - (1.96*lm2$std.error[4]))) 
com3 <- c("Whataboutism\nvs. Control", lm3$coefficients[4], (lm3$coefficients[4] + (1.65*lm3$std.error[4])), (lm3$coefficients[4] - (1.65*lm3$std.error[4])),
          (lm3$coefficients[4] + (1.96*lm3$std.error[4])), (lm3$coefficients[4] - (1.96*lm3$std.error[4]))) 
com <- data.frame(rbind(com1, com2, com3))
colnames(com) <- c("treatment", "effect", "cih90", "cil90", "cih95", "cil95")
level_order <- c("Identity Threat\nvs. Control", "Whataboutism\nvs. Identity Threat", "Whataboutism\nvs. Control") 

com$effect <- as.numeric(com$effect)
com$cih90 <- as.numeric(com$cih90)
com$cil90 <- as.numeric(com$cil90)
com$cih95 <- as.numeric(com$cih95)
com$cil95 <- as.numeric(com$cil95)
coma <- com


####results with controls####

lm1 <- lm_robust(list_items ~ as.factor(trt) + education + income + newjob + as.factor(trt)*sens, data=data1)
lm2 <- lm_robust(list_items ~ as.factor(trt) + lr_1 + natid_covar1_1 + income + newjob + as.factor(trt)*sens, data=data2)
lm3 <- lm_robust(list_items ~ as.factor(trt) + lr_1 + unemployed + as.factor(trt)*sens, data=data3)

com1 <- c("Identity Threat\nvs. Control", lm1$coefficients[6], (lm1$coefficients[6] + (1.65*lm1$std.error[6])), (lm1$coefficients[6] - (1.65*lm1$std.error[6])),
          (lm1$coefficients[6] + (1.96*lm1$std.error[6])), (lm1$coefficients[6] - (1.96*lm1$std.error[6]))) 
com2 <- c("Whataboutism\nvs. Identity Threat", lm2$coefficients[8], (lm2$coefficients[8] + (1.65*lm2$std.error[8])), (lm2$coefficients[8] - (1.65*lm2$std.error[8])),
          (lm2$coefficients[8] + (1.96*lm2$std.error[8])), (lm2$coefficients[8] - (1.96*lm2$std.error[8]))) 
com3 <- c("Whataboutism\nvs. Control", lm3$coefficients[6], (lm3$coefficients[6] + (1.65*lm3$std.error[6])), (lm3$coefficients[6] - (1.65*lm3$std.error[6])),
          (lm3$coefficients[6] + (1.96*lm3$std.error[6])), (lm3$coefficients[6] - (1.96*lm3$std.error[6]))) 
com <- data.frame(rbind(com1, com2, com3))
colnames(com) <- c("treatment", "effect", "cih90", "cil90", "cih95", "cil95")
level_order <- c("Identity Threat\nvs. Control", "Whataboutism\nvs. Identity Threat", "Whataboutism\nvs. Control") 

com$effect <- as.numeric(com$effect)
com$cih90 <- as.numeric(com$cih90)
com$cil90 <- as.numeric(com$cil90)
com$cih95 <- as.numeric(com$cih95)
com$cil95 <- as.numeric(com$cil95)
comb <- com

####combined figure####

coma$Controls <- "Without"
comb$Controls <- "With"
full <- rbind(coma, comb)

level_order <- c("Identity Threat\nvs. Control", "Whataboutism\nvs. Identity Threat", "Whataboutism\nvs. Control") 

ggplot(full, aes(x = factor(treatment, level = level_order), y = effect, color=Controls)) +
  geom_point(stat="identity", size=4, position = position_dodge(width = 0.15)) + 
  geom_linerange(aes(ymin = cil95, ymax = cih95), 
                 position = position_dodge(width = 0.15,),
                 size = 1) +
  geom_linerange(aes(ymin = cil90, ymax = cih90), 
                 position = position_dodge(width = 0.15),
                 size = 2) +
  ylab("Treatment Effect") +
  xlab("Treatment Type") +
  geom_hline(yintercept=0, linetype="dashed", color="blue") +
  theme_bw()+
  theme(
    axis.title=element_text(size=18), 
    axis.title.x = element_text(margin = unit(c(7, 0, 0, 0), "mm")),
    axis.title.y = element_text(margin = unit(c(0, 7, 0, 0), "mm")),
    plot.title = element_text(size=26),
    axis.text=element_text(size=16)
  )  

ggsave("figC1.png", width = 10, height = 6, dpi = 300)